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LONG-TERM GOALS 

The long-term goal of this work is to develop numerical simulation techniques for predicting wave- 
induced sediment transport and bed evolution in the coastal zone (order 10 meter depth outside the surf 
zone), at horizontal scales of tens of meters or less, and the scouring/burial around objects/obstacles 
(mine-like), partially buried in the bottom. This research was initially motivated by the Navy’s need to 
improve mine countermeasures (MBP), in support of joint littoral warfare, for which a key paradigm 
was to locate a clear or low mine- and obstacle-density path. More generally the laboratory validated 
modeling tools developed in this work will help improve both our fundamental understanding of, and 
ability to predict physical processes governing wave-current-induced morphodynamic changes, in 
coastal and estuarine environments of active sediment transport (e.g., ripple migration). 


OBJECTIVE 

In the past year, our objective was to extend the capabilities, integration, numerical efficiency and 
ultimately predictive abilities of two coupled numerical models simulating wave-(current)-induced 
flows and sediment transport over the seabed and around partially buried obstacles. These models 
were: (i) a two-dimensional (2D) fully nonlinear and inviscid Numerical Wave Tank (2D-NWT; 
initially developed at the University of Rhode Island; URI; Grilli and Subramanya, 1996; Grilli and 
Horrillo, 1997), and (ii) a three-dimensional (3D) fully viscous Navier-Stokes, Large Eddy Simulation 
model, with imbedded sediment transport model (3D-NS-LES; initially developed at Stanford 
University, and coupled in collaboration with Prof. Street and co-workers; Zedler and Street, 2001). In 
an earlier phase of this work, these models had been coupled and tested for a single idealized 
application, thus providing a proof of concept (Gilbert, 2005; Gilbert et al. 2005, 2007). The latter 
work clearly demonstrated the relevance of our modeling approach to simulate wave-induced flows 
and sediment transport around partially buried circular objects (mines) near the surfzone. Many 
aspects of the modeling, however, were initially oversimplified or idealized (i.e., flat bottom, quasi-2D 
object geometry, no moving bed algorithm, simple periodic incident waves) or arbitrary (sediment 
pick-up algorithm), and needed to be extended or refined for such computations to become both 
predictive and useful in solving scientific and practical problems of interest to this Navy program. Eor 
solving meaningful practical problems, the models also needed to be streamlined and optimized, and 
implemented and run on a massively parallel supercomputer cluster. 
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APPROACH 


In an earlier phase of this work, as part of collaborative work with Stanford University, the 2D-NWT 
developed at URI, based on Fully Nonlinear Potential Flow (FNPF) equations solved by a Boundary 
Element Method (BEM; Grilli and Subramanya, 1996; Grilli and Horillo, 1997), was extended and 
used to simulate the shoaling of incident periodic waves, including dynamic effects due to their local 
interaction with the bottom or an obstacle/mine, up to near breaking and subsequent dissipation in an 
absorbing beach (see sketch in Eig. 1). Such simulations were extensively validated by comparison 
with laboratory experiments (Grilli et ah, 1994, 1997a,b, 2004), which confirmed that ENPE theory is 
very accurate to model gravity wave transformation over complex bottom topography, up to and 
including wave overturning. 



Figure 1: Sketch of 2D-NWT, with piston wavemaker on the left-hand side, 
a sloping bottom, with a cylindrical, partially buried object, and 
an absorbing beach on the right hand side, where waves break. 



NWT Domain 


LES Domain 


Figure 2: Sketch of the NWT domain with embedded 3D NS-LES domain Qles 
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The hydrodynamic forcing from incident waves is calculated in the 2D-NWT, over a small near¬ 
bottom 3D region Qles containing the moving sediment covered seabed and the buried object, in 
which 3D-NS-LES and sediment transport computations are performed (Fig. 2). Assuming long- 
crested incident swells, wave forcing is expressed in terms of 2D velocity and pressure fields. In the 
earlier phase, computations in the 2D-NWT were made for a wide range of parameters (incident wave 
height and period, obstacle size and percentage of burial, and sandy bottom). Wave elevation and 
particle velocities near the obstacles were compared to laboratory experiments, performed at Arizona 
State University (Voropayev et ah, 2003), thus providing a physical validation of the simulated forcing 
wave fields (Grilli et ah, 2003, 2006). In the past year, we have implemented and tested the generation 
of more realistic fully nonlinear irregular waves with specified energy spectrum (see Results section). 

The NS-LES computations performed in Qles are fully 3D, which is essential for large-eddy 
simulations and sediment transport calculations. These computations are both initialized and forced in 
real time by wave-induced hydrodynamic fields. In the first year, only periodic swells were used and 
forcing for one period of NWT computations was repeated in a loop in the NS-LES model. In the past 
year, the model coupling was extended to apply to long time series of arbitrary wave forcing (see 
Results section). In the earlier phase, drawing on past experience in simulating oscillatory flows 
(Zedler and Street, 2001), after experimenting with several schemes, the 3D-NS-LES model was 
driven by specifying the wave-induced pressure gradients from the NWT model as source terms into 
the NS-LES momentum equations (Gilbert et ah, 2005, 2007). In the past year, a new formulation of 
model coupling has been implemented, based on a perturbation flow approach, which was shown to be 
both more consistent and accurate (see Results section). 

Because of the BEM solution, only the 2D-NWT boundary is discretized; the smaller Qles region (Fig. 
2), however, has a 3D finite difference numerical grid, highly-resolved in the vertical plane and less so 
in the direction normal to the main flow direction, to fully capture all the details of turbulent motion, 
particularly near the seabed, and accurately compute velocity, pressure, and Suspended Sediment 
Concentrations (SSC) (see, Zang et ah, 1993, 1994; Zedler and Street, 2001, 2002, for background of 
the 3D-NS-LES model). 

Rationale for this work ; The near-bottom wave-induced flow that forces sediment transport can be 
simulated using a standard 3D-NS-LES model, for arbitrary bottom topography with/without obstacles, 
at a few meter scale. This, however, requires computing at a larger scale (hundreds of meters) the 
background flow induced by nonlinear shoaling waves (and possibly currents), including effects 
induced by the bottom and obstacle/mine, from the far-field to the location of Qles- The 2D-NWT 
provides a very efficient and accurate solution to the larger scale problem and the full coupling NWT- 
LNS-LES, particularly using the most recent perturbation flow method, provides a seamless transition 
of computations from one model to another. The approach is also amenable to parallel computations. 

Key individuals at URI ; In the past fiscal year, the PI has worked on the project with two graduate 
students. Nate Greene, a new graduate student, joined the project in the Summer 2006 and has since 
concentrated on performing nonlinear irregular waves simulations in the NWT. Additionally, since the 
Spring 2006, Jeff Harris, a Ph.D student who has been awarded a 3-year NDSEG fellowship, has 
worked full time on all computational aspects of the project related to the 3D-NS-LES model, at no 
cost to ONR. 
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WORK COMPLETED 


In this two-year project, with the goal to perform more general and physically meaningful simulations 
of wave-induced sediment transport for seabed-obstacle evolution, we have pursued the development, 
testing and validation of the coupled NWT-NS-LES model along the following directions : 

(1) During the first year of the project, we had implemented and run the same 2D-NWT and 
3D-NS-LES models as developed in the previous phase, on UREs 8-node computer cluster. 
In the past year the model structure was modified in-depth by: (i) implementing a new 
perturbation flow approach for the NS-EES computations (further discussed in this report); 
(ii) and performing additional parallelizations and optimizations of both model codes. 

(2) Extension work initiated during the first year was continued in the past year, for the NS- 
EES and embedded sediment transport models, including implementing : (i) a sloping 
bottom; (ii) a more realistic and physically meaningful sediment pick-up algorithm (i.e., 
specifying the gradient of the SSC using a pickup function, instead of specifying the SSC at 
some point above the seabed); (iii) a logarithmic-layer near the bottom, rather than using a 
no-slip boundary condition (which is less accurate unless the viscous sublayer can be 
resolved); and (iv) forcing the model by arbitrary irregular, fully nonlinear, incident wave 
fields. Eor the latter extension (iv), limited simulations in the NWT had been made, but we 
now are able to generate fully nonlinear incident wave fields, corresponding to a specified 
energy spectrum (e.g., JONSWAP) (see Results section). 

(3) In the first year, the NS-EES code flowchart and structure had been modified to make them 
general in geometry and parameters, based on user-specified data (in the previous phase, 
the models were “hard-wired” to only solve one problem). In the past year, we run new 
applications aimed at testing and validating the various aspects of the model (see Results 
section in this report). More work is still required, however, to make the coupled model 
fully general in both incident wave fields and geometry of the obstacle and ocean bottom. 

(4) During the latter part of the past year, we were granted access to, and resources on one of 
the DoD massively parallel supercomputer clusters (with over 2000 processors). We have 
then been experimenting with both compiling and running the generalized and parallelized 
NWT and NS-EES models on this powerful computer. Once the coupled model is fully 
operational on this large cluster, we will run sensitivity analysis (i.e., convergence and 
accuracy) studies of model results to grid size and numerical parameter values. 

(5) During their initial development in the previous phase of this project, both models were run 
in sequence, which required a number of manual steps (some of these involving using small 
MATEAB programs). In the past year, we have achieve the fully automated coupling and 
running of both models on the clusters. 


RESULTS 

Wave generation : In the past year, we extended the 2D-NWT to simlate fully nonlinear irregular 
shoaling incident wave fields, with deep water wave input based on standard energy spectra (e.g., 
JONSWAP, PM). Eigure 3, for instance shows a 2D simulated free surface corresponding to a typical 
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New England shelf average seastate (Hs=1.2 m, Tp=5 s, depth 50 m). We are experimenting with two 
methods for performing such nonlinear irregular wave simulations in the NWT, with the coupled 
model, which are still being tested. In the first method, illustrated in Fig. 3, we : (i) initialize 
simulations by first creating a linear sea-state based on a specified spectrum, over constant depth, 
using the standard random phase method; (ii) use the efficient Higher Order Spectral (HOS) method 
(Dommermuth and Yue, 1987), with ramp-up procedure (Dommermuth, 2000), to evolve the linear 
sea-state towards full nonlinearity; and (iii) once reaching fully nonlinear quasi-permanent conditions 
in the HOS model over constant depth, waves are finally specified as input condition into the 2D- 
FNPF for further simulations. In the second method, illustrated in Fig. 4, we generate irregular waves 
in the NWT using a wavemaker, starting from still water conditions (e.g.. Fig. 1). The wavemaker 
stroke spectrum is obtained from the wave energy spectrum, using the usual linear superposition and 
transfer function. After generating waves, as in a physical wavetank, the wavemaker law of motion is 
iteratively corrected to achieve the target spectrum in the NWT. Although the second method is more 
computationally intensive, it allows for an arbitrary geometry of the NWT boundary, while the HOS- 
based method is limited to constant water depth. In Fig. 4, irregular waves are generated for a wind 
wave spectrum with Hs= 0.6 m, Tp= 2.2 s, in depth d = 6.4 m. We see a sequence of free surfaces 
generated in the NWT, with a ramp-up from Fig. 4a-d, to reach almost a steady state in Fig. 4e, where 
incident irregular waves are absorbed in the absorbing beach (for x ’> 4.5). The dominant wave height 
is around 0.\d, as expected and the dominant wavelength is about 0.5d = 3.2 m. 


J7(m) 



Figure 3: Sketch of fully nonlinear irregular waves computed with the HOS method, as input 
condition for the 2D-FNPF-NWT model, based on a JONSWAP spectrum 
(Hs=1.2 m, Tp=5 s, depth = 50 m). 


New 'perturbation flow’ coupling method : In this method, which we implemented in the past year, 
the total velocity and pressure fields (u, p) are divided into the sum of an inviscid wave part (Uw, Pv?} 
and a viscous perturbation part (Uv, pfj- The inviscid, irrotational, wave fields are known from the 
NWT solution, and can directly be computed in the BEM model for any time step, throughout the 3D- 
NS-FES domain Qles- NS-FES equations and boundary conditions are then rewritten in terms of the 
perturbation fields, while treating the wave fields as known forcing terms, at a given time step. The 
(linear) mass conservation equation keeps the same form as before for the perturbation velocity (since 
the wave velocity already satisfies continuity equation in the NWT), and is solved as a Poisson 
equation in the model. Because of nonlinear convective terms in the momentum equations, interaction 
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terms between wave and perturbation fields must be calculated in realtime. One advantage of the new 
formulation is, perturbation fields now naturally decay in the model, away from the seabed and 
obstacle, and naturally satisfy a radiation condition on the outer boundary of domain Qles- This was 
not the case in the earlier implementation (e.g., Gilbert et ah, 2007). Details of this formulation will be 
reported in Harris and Grilli (2007) and Grilli, Harris and Greene (2007). Note, this new 
implementation, with a flow decomposition, is exact in the equations and hence does not require any 
new approximation, except as before that the perturbation from the NS-LES flow on the NWT solution 
is neglected, which is relevant considering the scales of the wave and sediment transport flows. 
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Figure 4: Sketch of fully nonlinear irregular waves directly generated with a wavemaker (atx’ = 0) 
in the 2D-FNPF-NWT model, based on a JONSWAP spectrum (Hs= 0.6 m, Tp= 2.2 s, d = 6.4 m), 
with t/T = (a) 5, (b) 9, (c) 13, (d) 17, (e) 21, T/At = 100, and the absorbing beach starts atx’ = 4.5. 

We define x’ = x/d and z’=z/d. 
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Figure 5: wave-induced periodic flow (3 s period) around a 75% buried cylinder in a flat sandy 
bottom (with dso = 200pm). (a) Velocity vectors; and (b) sediment concentration. Axes are in meter. 
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Figure 6: Examples of theoretical flow computations for verification using the NS-LES model, (a) 
Blasius flow (middle section); (b) leading edge flow; (c) trailing edge flow. 


Flow and sediment transport simulations : In its current more general (i.e., using the perturbation 
flow) and parallelized version (see section on Work completed), the coupled NWT-NS-LES model can 
perform realistic and physically meaningful simulations of wave-induced sediment transport over a 
partially buried obstacle located on a sloping bottom covered with sediment. 

Last year’s report showed an example of eoupled NWT-NS-LES hydrodynamie model simulations for 
a purely periodic wave foreing. In the past year, we verified that these results could be retrieved using 
the more reeent model implementation, with the new eoupling method (flow perturbation) and parallel 
model implementation on the clusters. Ligure 5 shows an example of sueh verifieation computations 
for a wave-induced periodic flow around a 75% buried cylinder, in a flat sandy bottom (with dso = 200 
pm), at some time during one wave period. Velocity vectors are shown in Lig. 5a and sediment 
coneentration is shown at the same time in Lig. 5b. The foreing flow is a 3 s period wave, with 
maximum horizontal current velocity of 0.05 m/s. The arbitrary sediment pickup law is still used in 
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this computation, whereby the SSC is speeified at the nearest point above the seabed dependent on the 
frietion veloeity (see, Gilbert et ah, 2007). 

To better validate the various aspeets of the NS-LES model, in relation to numerical parameter values, 
partieularly regarding the flow indueed in the bottom boundary layer, in the past year, the model was 
run for standard or aeademie type flows, whose solution is either known from analytieal solutions, 
experiments, or other numerieal solutions. Figure 6 shows examples of such validation applications, 
for a Blasius flow (a) (i.e., viscous flow near the middle of a long, theoretieally infinite, plate, and 
flows at the leading (b) and trailing (e) edge of sueh a plate). The figure shows that eomputed flows, 
materialized by veloeity veetors, follow the expeeted pattern in the boundary layer near the plate. 

Other applieations are eurretnly being run to test the new sediment and moving bed algorithms; the 
model is thus expeeted to predict bed evolution, near and around a bottom obstacle, as a function of 
time (i.e., burial/seouring) for periodie or irregular wave foreing. Sueh results will be reported in the 
fortheoming papers (Harris and Grilli, 2007; and Grilli, Harris and Greene, 2007) 

IMPACT/APPLICATIONS 

We are developing an efficient and aeeurate numerieal simulation tool for large-seale far field wave 
foreing in a 2D-NWT model, and detailed representation of indueed near-bed flows and sediment 
transport, by large-eddy simulation, in a near field 3D-NS-LES model. This tool has great potential for 
both improving our fundamental understanding of wave-indueed sediment transport processes and 
predieting the burial and unburial of objeets, partally buried on the eoastal oeean bed, and for 
predicting ripple and dune formation and evolution (e.g., Voropayev et ah, 1999). 

In a larger eontext, this work will yield an experimentally-validated numerieal tool for studying 
nearshore wave elimates and indueed sediment transport and geomorphodynamie bottom ehanges. The 
effieieney and versatility of the proposed eoupled model will allow, in future work, investigating new 
elasses of problems of interest to various ONR programs. Speoifieally, the present problematie and 
methodology direetly relate to the objeetives of the ONR Wave/Mud Interaetion MURI and Tidal Flats 
DRI projeets. In particular, the model should ultimately help in assessing and simulating the effeets of 
dissipation, due to sandy or muddy bottom sediment and morphodynamie changes indueed on the 
bottom by the eombined effeets of waves and tidal currents, on nearshore wave transformations and 
wave parameters observable through remote sensing (e.g., wave steepness, wavelength, surfaee slope 
and eurvature). 
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